#include "../../headers/structsystem.h"
using namespace std;


int main()
{
    // command-line parameters
    cout << "Please enter P: ";
    double load_factor;
    cin >> load_factor;

    // create system
    StructSystem system = StructSystem(1);
    system.setJobname("BilinearTruss");

    // solve parameters
    double endTime = 100.0;
    double h = 1.0e-3;
    double zeta = 1;
    double length = 1;
    double theta = 45.0 / 180 * PI;

    // create particles
    Particle* p1 = new Particle(1, -length, 0);
    Particle* p2 = new Particle(2, 0, 0);
    Particle* p3 = new Particle(3, length, 0);
    Particle* p4 = new Particle(4, 0, length);
    system.addParticle(p1);
    system.addParticle(p2);
    system.addParticle(p3);
    system.addParticle(p4);

    // create material
    double E = 2.06e11;
    double yield_stress = 2.35e8;
    double Et = 2.06e10;
    double m = 0;
    UniBilinear mat = UniBilinear(1, E, Et, yield_stress, m);
    mat.setDensity(7850);

    // LinearElastic mat = LinearElastic(1);
    // mat.setE(E);
    // mat.setDensity(7850);


    // create section
    CustomSectionParas paras = {.A=1.0, .Sy=0, .Sz=0, .SHy=0, .SHz=0, .CGy=0,
                                .CGz=0.0, .Iyy=0.01, .Izz=0.01, .Iyz=0};
    CustomSection sect = CustomSection(1, paras);

    // create elements
    Link2DLD* e1 = new Link2DLD(1, p1, p4, &mat, &sect);
    Link2DLD* e2 = new Link2DLD(2, p2, p4, &mat, &sect);
    Link2DLD* e3 = new Link2DLD(3, p3, p4, &mat, &sect);
    system.addElement(e1);
    system.addElement(e2);
    system.addElement(e3);

    // constraints
    DOF c = {.key = {true, true, false, false, false, false},
              .val = {0, 0, 0, 0, 0, 0}};
    for (int i = 1; i <= 3; i++)
    {
        system.addConstraint(i, c);
    }
    system.info();

    // save model information
    string path = system.workdir() + "/" + system.jobname() + "/model";
    system.saveModel(path);

    StdArray6d f {};
    // f[1] = load_factor;
    system.addExternalForce(4, f);
    system.solve(h, zeta, true);
    system.setInternalForce();

    int nStep = ceil(endTime / h);
    int interval = ceil(1 / h);
    for (int j = 0; j <= nStep; j++)
    {
        system.solve(h, zeta, false);
        system.clearParticleForce();
        f[1] =  load_factor * (j * 1.0 / nStep);
        system.setExternalForce(4, f);
        system.setInternalForce();

        // save results
        if (j % interval == 0)
        {
            string path = system.workdir() + "/" + system.jobname() + "/" +
                          to_string(f[1]);
            system.saveParticleResult(path);
            system.saveElementResult(path);
            system.saveSupportReact(path);
        }
    }

    system.releaseContainers();
    return 0;
}